module mesh_manager

    use common_info
    use space_manager
    use time_manager
    use dynamic_quantities, only:       &
        folder_name => module_name,     &
        coordinate_name,                &
        quantities_name,                &
        num_quantity,                   &
        dynamic_quantities_is_initialized
    use dist_recon
    use file_io
    use message
    use run_manager

    integer :: num_grid(3)

    real, allocatable :: dx(:)
    real, allocatable :: x(:,:,:,:)
    real, allocatable :: f(:,:,:,:)
    real, allocatable :: dfdx(:,:,:,:)

    namelist /mesh_manager_parameters/ &
        num_grid

contains

    subroutine mesh_manager_init(namelist_file)
        character(*), intent(in) :: namelist_file

        integer :: i, j, k, l, o(3)

        ! module dependencies
        if(.not. dynamic_quantities_is_initialized()) then
            call message_show(ERROR_MSG,     &
                "mesh_manager",              &
                "mesh_manager_init",         &
                """mesh_manager"" needs "    &
                // """dynamic_quantities"" " &
                // "to be initialized first")
            call run_manager_end_run
        end if

        ! default values
        num_grid = 1
        num_grid(1:num_dim) = 50

        open(10, file=namelist_file)
        read(10, nml=mesh_manager_parameters, end=10)
10      close(10)

        allocate(dx(num_dim))
        allocate(x(num_grid(1), &
                   num_grid(2), &
                   num_grid(3), &
                   num_dim))
        allocate(f(num_grid(1), &
                   num_grid(2), &
                   num_grid(3), &
                   num_quantity))
        allocate(dfdx(num_grid(1), &
                      num_grid(2), &
                      num_grid(3), &
                      num_quantity))
        dx = (boundary(2,:)-boundary(1,:))/num_grid

        do k = 1, num_grid(3)
        do j = 1, num_grid(2)
        do i = 1, num_grid(1)
            o = (/i,j,k/)
            do l = 1, num_dim
                x(i,j,k,l) = boundary(1,l)+dx(l)/2+(o(l)-1)*dx(l)
            end do
        end do
        end do
        end do
        
        return
    end subroutine mesh_manager_init

    subroutine mesh_manager_parcel_to_grid
        integer :: i, j ,k, l
        character(20) :: label

        do l = 1, num_quantity
        do k = 1, num_grid(3)
        do j = 1, num_grid(2)
        do i = 1, num_grid(1)
            call dist_recon_evaluate(l, x(i,j,k,:), f(i,j,k,l))
        end do
        end do
        end do
        end do
        
        ! output
        if(time_manager_is_first_step()) then
            do i = 1, num_dim
                call file_io_create_data_card(output_file, &
                    folder_name, "mesh_"//trim(coordinate_name(i)))
                call file_io_output("constant", output_file,        &
                    folder_name, "mesh_"//trim(coordinate_name(i)), &
                    x(:,:,:,i), 3, num_grid)
            end do
            do i = 1, num_quantity
                call file_io_create_data_card(output_file, &
                    folder_name, "mesh_"//trim(quantities_name(i)))
            end do
        end if

        label = time_manager_cur_time_string()
        do i = 1, num_quantity
            call file_io_output(label, output_file, folder_name, &
                "mesh_"//trim(quantities_name(i)), f(:,:,:,i), 3, num_grid)
        end do

        return
    end subroutine mesh_manager_parcel_to_grid

end module mesh_manager
